Revealing endogenous conditions for Peto’s paradox via an ordinary differential equation model

Cancer, a disease intimately linked to cellular mutations, is commonly believed to exhibit a positive association with the cell count and lifespan of a species. Despite this assumption, the observed uniformity in cancer rates across species, referred to as the Peto’s paradox, presents a conundrum. Recognizing that tumour progression is not solely dependent on cancer cells but involves intricate interactions among various cell types, this study employed a Lotka-Volterra (LV) ordinary differential equation model to analyze the evolution of cancerous cells and the cancer incidence in an immune environment. As a result, this study uncovered the sufficient conditions underlying the absence of correlation in Peto’s paradox and provide insights into the reasons for the equitable distribution of cancer incidence across diverse species by applying nondimensionalization and drawing an analogy between the characteristic time interval for the variation of cell populations in the ODE model and that of cell cycles of a species. Supplementary Information The online version contains supplementary material available at 10.1007/s00285-024-02123-7.


Background on Peto's paradox
Cancer continues to be a perplexing issue within the realms of biology and medicine, prompting humanity's ongoing quest for understanding and resolution.As a disease intrinsically connected to somatic mutations, one might expect larger and longer-lived organisms to be more susceptible to cancer than their smaller and shorter-lived counterparts, given the increased number of cell divisions and the assumption that mutation rates remain constant across species during these divisions.However, this anticipated association is not observed, giving rise to the phenomenon known as "Peto's Paradox" (Peto et al. 1975;Peto 2016).This puzzling lack of correlation has been substantiated through a comprehensive investigation involving 191 species of zoo mammals (Vincze et al. 2021), further emphasizing the enigmatic nature of cancer incidence and its relationship with an organism's cell count and longevity.

Objective of the study
The intriguing question of how larger and longer-lived organisms have evolved mechanisms to counteract the increased likelihood of cancer development has piqued scientific interest, as unravelling this mystery may hold the key to advancing human cancer treatment.Prior investigations into Peto's paradox have explored factors such as tumour suppressor genes and energy turnover rates.For example, the comparatively lower cancer rates in elephants have been attributed to their additional copies of the TP53 gene (Abegglen et al. 2015;Sulak et al. 2016).Furthermore, Maciak and Michalak's research (2015) posits that larger animals' reduced energy turnover rates contribute to a decreased risk of cancer initiation.Salazar-Bañuelos' study (2019) examined the immune system, specifically the balance between reactive and suppressive lymphocyte clones, and attempted to elucidate the paradox using a Polya's Urn model.However, this approach yielded a negative correlation rather than the anticipated non-correlation between cancer risks and species' cell numbers.In our current research, we aim to shift the focus from clone ratios to the interplay within the cancer-healthy-immune system triad.By examining these interactions, we aspire to mathematically identify the conditions that give rise to the observed non-correlation between cancer incidence and species' size and lifespan, as the immune system's role extends beyond its internal workings.

Structure of the paper
This paper is organized as follows: Sect. 2 devotes to a literature review of the current research.Section 3 outlines the cancer-healthy-immune interaction model, while Sect. 4 delineates the requisite conditions for the non-correlation observed in Peto's paradox.Section 5 offers a visual representation of the findings and a discussion based on the result, and Sect.6 encompasses a conclusion of the findings, significance and future pathways.

Previous research on Peto's paradox and gaps
Peto's paradox refers to the observation that different species have a similar incidence of cancer despite their varying lifespans and cell numbers.There is a wealth of research aimed at identifying the underlying rationale that explains this non-correlation.Previous studies on Peto's paradox can be classified into three categories: 1. Specific genes: Some studies have focused on cancer-suppressing genes (Abegglen et al. 2015;Sulak et al. 2016), such as TP53, which is found in elephants and plays a crucial role in DNA repair and cell cycle regulation.However, this explanation may not account for such a universally existing phenomenon.2. Cell division and energy turnover rates (Maciak and Michalak 2015): Larger animals may have lower rates of cell division and energy turnover, resulting in a lower probability of developing cancer.Mutation rate considerations (Calabrese and Shibata 2010;Caulin et al. 2015;Nunney 2020) can also be classified into this category.3. Cell number considerations: These studies involve the "lethal tumour amount" (Nagy et al. 2007) and the stability or efficiency of the immune system (Salazar-Bañuelos 2019).The concept of an "immune escape threshold" (Palmer et al. 2018) has also been proposed in cancer research, though not specifically in the context of Peto's paradox.
We found that the previous research often neglected the relationship between species' lifespan and cell number, hence lacking a quantitative relationship between the influencing factors and lifespan or cell number.Our study seeks to fill these gaps by mathematically deriving the relationships between the influencing factors and lifespan or cell number.Recognizing that measurements on animals are difficult and data are scarce, our study tries to explain how these relationships can be derived through mathematical theory.
Several mathematical studies have examined Peto's paradox from various perspectives, including the energy supply and cancer waiting time viewpoint (Kempes et al. 2020), as well as the probability multistage model based on the driving mutation numbers (Calabrese and Shibata 2010;Nunney 2020).In contrast to these investigations, our focus will solely be on the cellular level, specifically addressing aspects such as self-growth and cell-cell interactions.

The role of the immune system in cancer development
Extensive research has been conducted on the role of the immune system in cancer development and treatment, making it a well-established concept.When analyzing cancer, the immune system cannot be disregarded.Since cancer development is not an isolated process within an organism, it is crucial to consider its interaction with other factors.Additionally, the aging of the immune system is widely recognized as one of the primary factors contributing to cancer (Palmer et al. 2018;Salazar-Bañuelos 2019).
In this article, we adopt the concept proposed by Palmer et al. (2018) that the effectiveness of the immune system acts as an upper threshold, referred to as the immune escape threshold (IET), for cancer development.We assume that this effectiveness, or competence, is proportional to the T-cell production and declines over time (Xie et al. 2017;Zhang et al. 2021).The reason we focus solely on T-cells lies in that CD8 + T cells are the most prominent anti-tumour cells and cause the direct destruction of target cells (González et al. 2018;Whiteside 2022), and tumour infiltrating cells, which are immune cells that have migrated from the blood stream into tumour tissues with the purpose of tumour elimination, consist mostly of CD3 + T cells and a smaller proportion of B and NK cells, whose functions are less clear (Liu et al. 2023).Thus, the immune system is associated with the cell count of the organism and can be treated as a function of time.Further details of the model will be presented in Sect.3.

The modelling of cancer growth and immune response
In the realm of mathematical biology, previous studies have utilized ordinary differential equation (ODE) systems to model the interactions among cancerous, healthy, and immune cells, aiming to comprehend cancer growth and the immune response (Alharbi and Rambely 2020;Dritschel et al. 2018).Rather than solely examining the specific changes in particular types of immune cells, we emphasize the characterization of the overall effectiveness of the immune system in combating cancer.

Theoretical framework
Cancer development within an organism extends beyond the mere presence of cancer cells; it is a multifaceted outcome resulting from the interplay between cancerous and non-cancerous cells.To adequately capture the overall cancer progression, we employ a Lotka-Volterra (LV) ordinary differential equation (ODE) model.
Originally used to model the change of prey and predator populations over time, the LV model is well-suited for capturing the dynamics of interacting populations over time.Generally, the model has the following form: Here, r 1 represents the natural growth rate of the prey species, and r 2 represents the death rate of the predator species in the absence of predation.Coefficient α is the decay rate of the prey species due to predation, while β is the growth rate of the predator species.The right-hand side of the LV equation can be divided into two components.The first component describes the self-growth and death rate, while the second component characterizes the gain or loss resulting from interactions with other species.This structural framework allows for the extension of the LV model to a broader range of applications, accommodating a greater number of species with their self-changes and interactions, including competition, cooperation, and more.The process of cancer progression, involving cell growth, death, and competition, is a suitable candidate for the LV model.

Model components
Striving for a simple yet comprehensive representation of cellular activity, we incorporate the following key players as shown in Fig. 1: 1. Cancer cells, which exhibit natural growth constrained by their carrying capacity, vie with healthy cells for resources and face immune cell attacks.2. Healthy cells, as the precursors of cancer cells, experience natural growth under their carrying capacity while competing with cancer cells for resources.3. Immune cells whose function is subject to aging and both activation and impairment caused by interactions with cancer cells.
Previous research has explored the aging of the immune system as a fundamental factor contributing to age-related diseases (Palmer et al. 2018;Salazar-Bañuelos 2019).In our model, the initial cellular values determine the immune system's efficacy, which, in turn, sets the upper limit for controlling cancer cell proliferation.If cancer cell numbers exceed this threshold, we assume that immune escape is to occur (Palmer et al. 2018).
Consider the following ODE system: Fig. 1 The relationships between the model components, namely cancer cells, healthy cells and the immune system here, C, H and I represent the populations of cancer and healthy cells, and the immune threshold respectively.Coefficients r 1 , r 2 and r 3 represent the natural growth rate of cancer cells, healthy cells, and the age-related decay of immune cells.We assume that the immune ability is linearly dependent on the number of immune cells so that the aging of the immune system would decrease this immune escape threshold monotonically with age (Palmer et al. 2018).Carrying capacities of cancer and healthy cells are denoted as K 1 and K 2 respectively, while Greek parameters α, β, γ , δ are the coefficients for cell interactions.In particular, coefficient α includes the overall effects of competition with and generation from healthy cells to cancer cells.We assume that the effect of the competition is stronger than that of canceration so that the net contribution to cancer cells is negative.Coefficient β corresponds to the loss of cancerous cells by the killing of immune cells, coefficient γ describes the strength of healthy cells' competition with cancer cells, and coefficient δ controls the net consumption of immune capacity when fighting cancer cells.It represents the overall effect of activation and depletion by cancer cells.We assume δ > 0 to represent a weakened health condition during the disease (Zheng et al. 2021).There is no direct interaction between healthy and immune cells here.Initial amounts of the three kinds of cells are denoted as C 0 , H 0 and I 0 .Because d I dt is only influenced by negative terms, I 0 should be the highest immune threshold, which corresponds to an intact immune system.

Cancer incidence
Let us assume that all species would have such an ODE system but are distinct by their different lifespanT 's.For example, for a speciesn, the dynamics of cell numbers can be described by the LV model with a unique set of parameters, namely,r 1 (n) , r 2 (n) , α (n) , β (n) , ..., T (n) , H 0 (n) , I 0 (n) .The initial healthy and immune cell amounts H 0 (n) , I 0 (n) correspond to the species' cell number mentioned in Peto's paradox.On the other hand, the population's average immune escape time is represented byt c (n) , whereC (n)  (n) , the count of cancer cells will be always smaller than the immune threshold, therefore we assume a null incidence of cancer.For timet > t c (n) , the cancer cells will surpass the immune threshold, implying a higher risk for tumorigenesis.During this period, we assume that there is a constant probability q for occurrence of cancer within a unit time increment.Therefore, the cancer incidence throughout the lifespan can be represented by 1 Apparently, a larger (smaller) value of t c (n) T (n) indicates a delayed (accelerated) emergence of unrestrained growth of cancer cells throughout the lifespan, which relates to a lower (higher) incidence of cancer.If values of t c (n) T (n) are the same for all species, these species would exhibit identical cancer incidence, which aligns with the observation of Peto's paradox.More specifically, for the non-correlation between cancer incidences and species to hold, the crossing points of C (n) (t) and I (n) (t) are expected to have the following relation: Note that prior studies have demonstrated a power law relationship between species mass and lifespan (Speakman 2005;Vazquez 2017), though previous research on Peto's paradox often neglected the correlation between species' cell numbers and lifespans.Assuming that each cell in different species has similar masses and that differences in the body mass mainly arise from cell numbers rather than cell sizes (Savage et al. 2007;Schmidt-Nielsen 1984), we can postulate a power law relationship between species' cell numbers and their lifespans.In particular, we assume the following relations for a specific species, where a, b, c and m are constants.Furthermore, we assume that initial cell numbers are proportional to the carrying capacities, namely.
1, which indicates a small number of cancer cells and almost a full amount of healthy cells at the initiation.Accordingly, initial cell numbers are quantitatively connected to T .

Stability properties of the ODE model
Stability analysis shows that there are four equilibrium points (C * , H * , I * ) for this LV system, i.e.

Equilibrium point
To ensure the existence of the fourth equilibrium point and keep the cell numbers lower than the carrying capacities (0 At the limit of these conditions, we have the following properties: 0, the fourth point becomes the third point (K 1 , 0,0).
points which form a line as: The stabilities for these points are summarized as follows, while the proof can be found in Supplementary Material S2.
1.The first point (0,0, 0) is always unstable.2. When r 1 − α K 2 < 0, among the first three equilibrium points, only the second point (0, K 2 , 0) is stable.3. When r 2 − γ K 1 < 0, among the first three equilibrium points, only the third point ) is always stable if the conditions for its existence are satisfied.

The tipping time
Since we represent the cancer incidence by the value of q 1 − t c T , at least we need to confirm the existence of this tipping time t c , or equivalently, the increasing C(t) and the decreasing I (t) should have an intersection point at a certain time within the lifespan of the species.
Settings of the parameters α, β, γ , δ, a(or K 1 ), b(or K 2 ), c(or I 0 ), r 1 , r 2 , r 3 , which involves some arbitrariness, will affect the existence of the intersection point.Although a thorough investigation of the value range for all the parameters may not be feasible, we may set some requirements for setting these parameters to recover those essential biological properties.These requirements include the following items: 1.The exponent of the power law is set as m 4.76, according to Speakman's study (2005).2. The tipping time t c should be found more likely in the late stage of life.Suppose that q takes the value of 0.5, for human cancer incidence t c /T would be around 0.992.For a 0.2 cancer incidence rate, t c /T equals 0.6.For easiness of comparison later, we let t c T > 0.6.3.There should be no significant decrease of healthy cells before the tipping time, hence H (t) const for t < t c .Here we allow a 5% decrease of healthy cells from H 0 .4. The effective depletion rate of the immune cells by cancer cells should be larger than the decay rate, i.e., K 1 δ > r 3 .
Our strategy is to let requirements 1, 4 be satisfied firstly, then we further tune the parameters manually to fulfill requirements 2 and 3. To show that we can find these parameter settings as well as the reasonableness of these requirements, we give five visual examples for the comparison of the effects of parameter settings.
• Case I: tipping time exists and both requirements 2 and 3 are met.The basic parameters are set as T 1, b 10000, a c 5000.Parameter δ is set to 0.001 to make requirement 4 satisfied, i.e., K 1 δ > r 3 where r 3 2.We adjust parameters α, β, γ so as to make their values close to each other, resulting in α 0.0001, β 0.00015, γ 0.00012.This will also ensure that r 1 > αK 2 , r 2 > γ K 1 hold for r 2 2, r 1 10.Under such a parameter setting, we may see in Fig. 3 that the intersection locates around t c 0.875, and the healthy cells barely change before t 0.85.• Case II: tipping time exists but both requirements 2 and 3 are not met.In this case, we set the parameters as the following: r 2 1000, r 1 5r 2 5000, α 0.35, β 0.2, γ 0.15.We can see from Fig. 4 that the intersection point locates around t c 0.22, in the meantime, both the healthy and immune cells start to decrease from the beginning of the simulation.
• Case III: tipping time exists but requirement 2 is violated.In this case, we set the parameters as the following: r 2 10, r 1 5r 2 50, α 0.0005, β 0.0006, γ 0.0002 which also let r 1 > αK 2 , r 2 > γ K 1 be satisfied.We can see from Fig. 5 that the intersection point locates around t c 0.22, which is quite early, and the large change in the count of healthy cells occurs around t 0.2, after which it becomes almost constant.
• Case IV : tipping time exists but requirement 3 is violated.In this case, we set the parameters as the following: r 2 1000, r 1 5r 2 5000, α 0.46, β 0.11, γ 0.11.We can see from Fig. 6 that the intersection point locates around t c 0.88, however, the count of healthy cells starts a slow but persistent decrease from t 0.37.• Case V : tipping time does not exist.In this case, we set the parameters as the following: r 2 1000, r 1 5r 2 5000, α 0.46, β 0.2, γ 0.15.We can see from Fig. 7 that the intersection point of C(t) and I (t) does not exist throughout the lifespan.
We see that in case V, β and γ are taken the same value as those in case II and that a larger α will delay the intersection, as indicated in the ODE equations.By comparing case IV and case V, it is obvious that a smaller β would bring ahead the intersection.In Case III and IV, we are also able to see the sudden growth in cancer cell count within a small period of time but not from the beginning of the simulation, representing the fast progression of tumour and the small cancer incidences at young ages.

Sufficient conditions for Peto's paradox
Recall the model from Sect.3.2, and suppose that for each species n there would be a system of ODEs, dC (n)  dt r We can represent the cell growth rates and interacting coefficients as equations of T (n) , so that they are in relation with the lifespans: where d (n) , e (n) , f (n) , g (n) , h (n) , k (n) , l (n) are positive constants.The d (n) , e (n) , f (n)  are unitless and g (n) , h (n) , k (n) , l (n)  First we would like to nondimensionalize the systems.Defining the characteristic scales of cell numbers and time as follows.
The dimensionless cell numbers and their initial values can be defined as the following using the characteristic scales.
Substituting Eq. ( 10) into the ODE system, we can derive the dimensionless system by rescaling the time as τ t t s , which gives where the dimensionless coefficients are Substituting the characteristic quantities and Eqs.(7-9) into dimensionless coefficients gives Meanwhile, the dimensionless intersection time is so the cancer incidence after nondimensionalization can be represented as We then propose the sufficient condition for Peto's paradox.

Proposition:
The following relations are a sufficient condition for Peto's paradox to occur.d (n) d, e (n) e, f (n) f, g (n) g, h (n) h, k (n) k, l (n) l where d, e, f, g, h, k, l are constants, for each species n.
Next, we would like to prove this proposition.
Proof When d (n) d, e (n) e, f (n) f, g (n) g, h (n) h, k (n) k, l (n) l for each species n, the dimensionless coefficients become constants, i.e., f .(15*) Then the dimensionless system can be rewritten as Obviously, the intersection time t c obtained from Eqs. (11*-13*)become a universal quantity across these species, i.e., Thus, Eq. 6 is satisfied and each species would have the same cancer incidence of q 1 − t c T .
Note the relationships in the Proposition can be rewritten as The presence of a common cancer incidence in the dimensionless ODE system indicates the emergence of Peto's paradox under the condition that biological properties such as carrying capacities, strengths of cell interactions and growth rates are solely determined by the expected lifespan of the species.Here we performed a series of numerical experiments to show the emergence of Peto's paradox in the ODE model.In the experiment, for any species with an expected lifespan T , we first identify the tipping time t c at which the population of cancerous cells surpasses that of immune cells, afterwards we further show that the representative cancer incidence q(1 − t c T ) s are the same for different species.
Let us assume three species have expected lifespan as T 1,10,100, with the evolution of cells falling under Case I of Sect.4.2.The interacting coefficients and growth rates are set following Eqs.(7, 8*, 9* ), where e 2, f 2, d 5e, b 100000, a 0.5b, c 0.5b, g 0.00001, h 0.000015, k 0.000012, l 0.0001.We use Matlab and the 'ode45' function for numerical simulations with default relative error tolerance(RelTol) and absolute error tolerance(AbsTol) and 'NormControl' 'on' in 'odeset' for the ode function setting options.Figure 8 illustrates the time evolution of cell numbers over the lifespan of the three species.Clearly, there exists a tipping time at which the decreasing population of immune cells has a "dead cross" under the increasing cancer cell population for each species.Now if we rescale the time by species' lifespan, and also rescale cancerous and healthy cell populations by their carrying capacities, and immune cell population by their initial value, we can see in Fig. 9, that the dimensionless tipping times, or equivalently the representative cancer incidences are the same for the three species.
Through the nondimensionalization and the above numerical experiment, we may conclude that Peto's paradox emerges in the proposed ODE system if the relationships between species' lifespans and their biological properties including numbers, interaction strengths and growth rates of cells satisfy a set of power laws shown in Eqs.(7, 8*, 9*).These relationships are the conditions for the noncorrelation of cancer incidence with the body of mass stated in Peto's paradox, identified from the ecological interactions of cancer, healthy and immune cells.
As shown in Eqs.(14* and 15*), growth rates and interacting coefficients become constants in the scaling of equations.Empirically, as long as the products of scaling constants, namely g × b, h × c, k × a and l × a are kept as constants, whatever the specific values that a, b, c, g, h, k, l took, the intersection point would always stay the same.For a comparison of this numerical experiment and Case I in Sect.4.2, a, b, c are set 10 times larger, and g, h, k, l are taken 10 times smaller, while the intersection time remains unchanged.This fact may provide us with a high degree of freedom and easiness if we calibrate the LV system with a specific species in the future when biological data of cell growth and cell interactions are available for that species.

Sensitivity analysis on the initial condition
As explained in Sect.3.3, initially cancer cells should be the minority.Hence we set C 0 K 1 , where is a small number.In the numerical experiments to find the cancer incidence, typically we set 0.0001.Certainly, this setting will also affect the intersection point of cancer and immune cell populations.We conducted a sensitivity analysis on this parameter by varying the value within ±2% in Experiment I and set a large value 0.001 for this parameter.Other parameters, r 1 , r 2 , r 3 , α, β, γ , δ are unchanged to satisfy the conditions for the intersection.The results of this sensitivity analysis are shown in Table 1.
From these experiments, we find that a larger value of will result in an earlier intersection of cancer cells and the immune threshold, therefore a higher cancer incidence.Nevertheless, since the parameters satisfy the noncorrelation conditions, results will not vary for all the species satisfying the scaling relations between their biological parameters and their expected lifespans.

A comparison with previous research
To better understand the findings from our numerical experiments, we give a brief comparison of our results with previous studies before we provide a deeper discussion in Sect. 5. Previous research suggested that the slower division rates of larger animals could account for Peto's paradox (Maciak & Michalak 2015).In our research, the sufficient condition for the non-correlation of cancer incidence across species says that larger species may have slower natural growth rates, which would be inversely proportional to the species' lifespan or inversely proportional to m √ S if S represents the size measured by cell numbers.In the case of the immune system hypothesis, the age-related decay of the immune system is linearly proportional to the lifespan in the sufficient condition.
Additionally, our model examines cell interactions beyond the isolated growth of different types of cells.What we found is that the changing rates due to interaction are inversely proportional to S × T or (∼ T m+1 ), the product of the size and lifespan, as part of the sufficient condition.Unlike the approaches adopted by multistage models (Calabrese & Shibata 2010;Nunney 2020), our model represents cancer incidence through looking at the intersection time of the number of cancer cells and the immune escape threshold.

Interpretation of the non-correlation conditions
Recall that when species satisfy the conditions stated in Eqs.(7, 8*, 9*), the nondimensional tipping times for all species would be the same, making the non-correlation between the incidence and species.These equations mean that for different cells in all the species, the growing coefficients are proportional to T −1 , and the interaction coefficients are proportional to T −(m+1) .Notice that T m+1 represents product of lifespan and size of cell numbers.

Understanding the non-correlation conditions
Let us assume that there are three species A,B, and C differentiated by their longevity, namely T 1 : T 2 : T 3 1 : 2 : 3.Moreover, the power law relation between cell numbers and lifespan (Eq.( 7)) holds for all these species.If we consider the simplest case of m 1, the ratio of body size is also S 1 : S 2 : S 3 1 : 2 : 3.In Fig. 10, we draw three cuboids to represent these species.In particular, the width of a cuboid represents its body size, and the length of a cuboid represents its lifespan.The small rectangular region in a cuboid denotes a unit time interval for cell behaviours, e.g., growth and interactions.
According to the discussion in Sect.5.1, the growth rates of cancer, healthy and immune cells of a species are inversely proportional to its lifespan, and the interaction rates are inversely proportional to the product of lifespan and the number of cells.Fig. 10 An illustration for the inverse proportional growing rates to the lifespan and inverse proportional interaction rates to the multiplicity of the lifespan and the size of cell numbers, but the same number of "intervals" with species A,B,C from top to bottom Thus, we may infer that the ratios of time intervals for different species must also follow a 1:2:3 pattern (i.e., the ratio of lifespans) so that all species must have the same number of time intervals in their lifespan bars.
Mathematically, the change of cell number of type i between the k th and (k + 1) th interval for a particular species whose lifespan is T can be calculated as the following where t T is the interval length for species with lifespan T , c 1 T N i T (t) characterizes the self-growth of the cells, and c 2 T m+1 N j T (t) characterizes the interaction with another cell type, and c 1 , c 2 are constants.For three species with lifespan ratio T 1 : T m according to Eq. ( 7).With N T i, j (k) proportional to T m , the ratio of the RHS of Eq. ( 16) is have ratio T 1 : T 2 : T 3 .Therefore, all species have the same number of intervals, or the number of cell cycles for different species are the same, if their body sizes are proportional to their lifespans (following Eq. ( 7)), and their coefficients of growth and interaction depend on their lifespans accordingly (following Eqs.(8*, 9*)).

Potential mechanism and related evidence
Notice that we have tentatively concluded that each species possesses an equal number of intervals.However, the precise definition of an interval in real-life remains ambiguous, although it presents an intriguing issue for investigation.Specifically, we are interested in determining which factors, such as growth and interaction rates, lifespan, and the number of intervals, influence one another.In terms of the self-growing factors, the number of intervals may be linked to Hayflick's limit.Nevertheless, it remains unknown whether different species share a common division limit.Regarding cell interactions, it is uncertain whether different species exhibit a universal frequency of interactions.Meanwhile, the length of telomeres, which is associated with Hayflick's limit and aging, could be a potential factor contributing to the absence of a correlation between cancer incidence and lifespan.It is widely recognized that the rate of telomere shortening serves as a predictor for lifespan (Whittemore et al. 2019).Comparable patterns can be observed among species in the biological world, albeit not necessarily at the cellular level.For example, despite variations in the lifespan of mammals, they tend to have a similar number of heartbeats throughout their lives (Levine 1997).
Furthermore, the study conducted by Cagan et al. (2022) demonstrated an inverse relationship between the mutation rates per genome per year and lifespans, with the mutation burden (total number of mutations per genome) being relatively consistent across different species.
As demonstrated in Eqs.(8* and 9*), maintaining a constant product of growth rates and lifespans and a constant product of interaction rates, lifespans, and cell numbers is sufficient for the emergence of Peto's paradox.However, obtaining direct evidence from real-life scenarios to determine whether such a constant value exists poses a challenge.We can only provide estimates for some of these constants based on the limited data available.The growth rate of cancer cells in a person is inferred from the doubling time of HeLa cells, which ranges from 17.5 to 32.3 h (Tang 2019).Similarly, the growth rate of healthy cells in humans (Homo sapiens) and mice (Mus musculus) is estimated using their respective doubling times.According to BioNumbers (Milo et al. 2010), the doubling time of human cells ranges from 24 to 72 h (search terms: 'doubling time of human cells ', BNID 103762, 112761, 102096, 109938), while the doubling time of mouse cells varies from 4 to approximately 20 h [BNID 109063, 110686, 113093].The growth rates of cells of several other species are also calculated by the doubling time of their finite cell lines, from the database Cellosaurus (Bairoch 2018) (search terms: 'species' name', 'doubling time', 'finite').The used data can be found in Supplementary Material S1.Moreover, Palmer et al. (2018) reported that r 3 of human is approximately 0.044 year −1 .Utilizing these values and the maximum lifespan values from the database AnAge (Tȃcutu et al. 2017) (these values are also recorded in Supplementary Material S1), the plausible values for constants in Eq. (8*) are estimated as follows:e ∼ 982; d ∼ 29872; f ∼ 5.39.
Nevertheless, it is evident that obtaining more precise estimates, validating these constants, and estimating others will necessitate more extensive datasets.

Other potential parameter sets
In Sect.4.3, the parameters we used (Case I) show a slow increase of cancer cells and a negligible decrease of immune cells.But certainly it is not the only parameter set.For another case, where b 10000, a 0.5b, c 0.085b, r 2 1000, r 1 5r 2 , r 3 2, α 0.476, β 0.6, γ 0.18, and δ 0.001, as shown in Fig. 11, we see a rapid increase in cancer cells and a sharp decrease in healthy cells within a short period.
The abrupt cancer progression in this case can be viewed as a first order phase transition that takes place with constant normal cell interactions and the aging of the immune system, without external environmental effects.Such characteristics exist in the parameter range β ∈ (0.462, 0.638).The study from Yu et al. (2015) reported the fast progression of pancreatic ductal adenocarcinoma from low-to-high tumour stages.

Limitations of our study
The interaction coefficient δ is positive in our study to represent a weaker health condition during the disease (Zheng et al. 2021).Future study may need to investigate the case where δ < 0 if necessary.In Eq. ( 7) we assumed the fraction of initial immune cells is the same across all species.Studies by Ruhs et al. (2020) and Downs et al. (2020) show the scaling relationships between body mass and lymphocytes in birds as log 10 lymphocyte concentration cells L 0.58 ± 0.17 + 0.003 ± 0.03log 10 (body mass(g)) and for 250 + mammalian species as log 10 lymphocyte concentration cells L 0.76 (0.51 : 0.92) − 0.04(−0.08: −0.02)log 10 (body mass(g)).
And according to Lindstedt and Schaeffer (2002), the relationship between animal blood volume and body mass is Blood V olume (ml) 71.5body mass(kg) 1.01 .
If applied to mammalian species, we have lymphocyte number u * body mass 0.97 u * T 4.6172 , where u is a constant.Then larger or longer-lived species might deviate more from this assumption (Eq.( 7)).However, the blood volume equation is obtained from only 4 species, necessitating further studies across different species.Moreover, the lymphocyte concentrations are not necessarily from very young ages.According to Bjornson-Hooper et al. (2022) and Jiao et al. (2024), the frequency of different immune cell types varies in humans, mice, and non-human primates or vertebrate species, yet whether part of the data are obtained from the very young age are not clear.Therefore, we are aware that more data are required for the justification of this assumption.
The estimation of d, e, f values in Sect.5.3 are based on the limited data available, so the accuracy can not be guaranteed.More accurate estimates should be made when more data become available.Additionally, whether the inverse relationship between growth rates and lifespans exists in real life is a good question to be asked.Current equations for fitting data on healthy cell growth rates versus lifespans are not necessarily inverse proportional functions, possibly due to low data density.Furthermore, the growth of cell lines may not accurately represent in vivo growth.
Interaction coefficients α, β, γ , δ in Eq. (9*) should be understood as the averaged rates over the entire population of a particular species.For example, the rate β for the interaction between cancer and immune cells can be calculated as the following (Kempes et al. 2020), where T i is the lifespan of a specific individual inside the species and β i (t) is the interaction coefficient of this individual at time t.Then the average is taken over all N individuals to give the species' average coefficient β.Therefore an ODE system with an average lifespan T represents a species in an average sense, and the non-correlation conditions for Peto's paradox are also to be interpreted statistically.Naturally, breakdowns of the non-correlation condition within a species can be readily observed due to factors such as gender, specific weights of individuals, and so on.Additionally, we cannot dismiss the possibility of individuals within a species experiencing significant fluctuations or perturbations.For example, in cases where new cancer cells emerge, the involvement of more immune cells leads to continued feedback resulting from this interaction.Moreover, we acknowledge that a larger immune system may possess more comprehensive and robust mechanisms.If considering that larger or longer-lived species might deviate more from our assumption of immune cell count (I 0 cT m ) and assuming the same proportions of different immune cell types across species, it is possible that these species have a more effective immune system.However, while the ODE model may not capture unique behaviours at specific periods or individual levels, we believe it is a suitable choice for studying different species with varying parameters.To further improve the ODE model, we may consider the fact that cancer cells C 0 starts from a small number instead of null initialization.If requiring C 0 0, we may modify the equation for cancer cells into where p represents the probability of mutational transition from a healthy to a cancerous cell and α here would only represent the effect of competition.Now at initiation, there would only exist the non-mutated type (healthy cells) and the environment (immune threshold).Note that p could be a function of p(t) over time, considering that cancer probability increases with age (Palmer et al. 2018).
Another main limitation of the ODE model is its omission of spatial changes as well as tumour migrations.To address this, a partial differential equation (PDE) model or an agent-based model (ABM) may provide a more accurate characterization of spatial behaviours.
For the use of the Lotka-Volterra model, the shortcomings lie in that one assumes a simplified interaction between populations, which may not fully capture the complexity of the tumour microenvironment.For example, the signaling pathways and feedback mechanisms are not accounted in the model.
Lastly, it is not our intention to solely attribute cancer development to intrinsic factors while disregarding other potential contributing factors.We are also aware of arguments suggesting the involvement of cell size (Maciak 2022;Maciak & Michalak 2015;Savage et al. 2007), and the need for additional considerations when establishing a linear relationship between cell number and animal size.Overall, we recognize that Peto's paradox is a complex phenomenon, and a combination of undiscovered factors may have contributed to this intriguing phenomenon.

Conclusion
In our study, we employed a Lotka-Volterra ODE model to identify the non-correlation conditions of Peto's paradox.Our findings indicate that the rates of cancer and healthy cell growth, as well as the age-related decay rate of the immune system, per unit time, should exhibit an inverse relationship with the species' lifespan for the sufficient condition.Furthermore, the interaction rate between cells per unit of time per influencing cell should be inversely proportional to the product of the cell number and lifespan.
While previous studies have explored the potential role of cell division rates and the significance of the immune system and its aging in Peto's paradox, few have integrated these factors to provide a comprehensive mathematical explanation that elucidates the specific relationship between these rates and lifespan.To our best knowledge, this paper is the first to present an ODE model between cells for understanding Peto's paradox and utilizes nondimensionalization to facilitate comparative oncology.
As discussed in Sect.5.5, future research should explore the spatial aspects of tumour migration and reduce dependence on initial conditions.Additionally, further investigation is required to determine whether and why each species exhibits similar division limits and the number of cell interaction times, and how this relates to other genetic and organismal levels.By addressing these aspects, we can gain insights into the equitable nature of cancer across various species.

Fig. 2
Fig.2The phase portrait near the C − H plane in the special case r 1 α K 2 , r 2 γ K 1

Fig. 3 Fig. 4 Fig. 5 Fig. 6 Fig. 7
Fig.3Case I.The upper panel is the 3D phase portrait of the system and the lower panel shows the dynamics of the count of cancer, healthy and immune cells.The trajectory of the time evolution of the system state is shown by colored symbols, where the color (from blue to red) is scaled with the time (from 0 to T ).Graphs in the lower part shows the intersection point of cancer and immune cells when b 10000, a 0.5b, c 0.5b, r 2 2, r 1 5r 2 , α 0.0001, β 0.00015, γ 0.00012 and δ 0.001

Fig. 8 Fig. 9
Fig. 8 Intersection point for C(t) and I (t) for T 1,10,100.The vertical axes represent cell numbers and the horizontal axes represent the lifespans.Here parameters of ODE model are set as the following: e 2, f 2, d 5e, b 100000, a 0.5b, c 0.5b, g 0.00001, h 0.000015, k 0.000012, l 0.0001

Table 1
Value of t c /T for